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SUMMARY 


A procedure is described for the automated optimum design of airplane wing 
structures subjected to multiple behavior constraints. The structural mass of the 
wing is considered the objective function. The maximum stress, wing tip deflection, 
root angle of attack, and flutter velocity during the pull-up maneuver (static load) , 
the natural frequencies of the wing structure, and the stresses induced in the wing 
structure due to landing and gust loads are suitably constrained. Both deterministic 
and probabilistic approaches are used for finding the stresses induced in the air- 
plane wing structure due to landing and gust loads. 

The procedure is illustrated with two examples. The first is a wing represented 
by a uniform beam with a cross section in the form of a hollow symmetric double 
wedge. The airfoil thickness and chord length are the design variables, and a graph- 
ical procedure is used to find the optimum solutions. The second example is a super- 
sonic wing represented by finite elements. The thicknesses of the skin and the web 
and the cross-sectional areas of the flanges are the design variables, and nonlinear 
programming techniques are used to find the optimum solution. 

For the examples considered, the optimum value of the objective function is 
higher when the landing and gust stresses are treated as random variables than when 
they are considered deterministic quantities. 


INTRODUCTION 

Considerable progress has been made in airplane structural optimization over the 
last two decades. Several investigators have considered both strength and aeroelas- 
tic requirements in airplane optimum structural design problems. In reference 1, a 
double-wedge airfoil is idealized as a bar, and a steepest descent method is used to 
minimize the total energy required to drive the airfoil through a specified flight 
condition. The behavior constraints considered in the optimization are the root 
angle of attack, the elastic deflection at the leading edge of the wing tip, the 
principal stress in the skin at the root of the wing, and the Mach number of bending- 
torsion flutter. The automated preliminary design of wing structures was considered 
in references 2 and 3 by idealizing the wing as a plate. In reference 2, the wing is 
idealized as an isotropic sandwich plate with a variable thickness cover, and an 
interior penalty function method is used to minimize the weight of the cover by vary- 
ing its thickness distribution. The method of feasible directions is used in refer- 
ence 3 to determine a simplified configuration of the wing structure, wing depth 
distribution, and cover panel thickness distribution by minimizing a linear combina- 
tion of the aerodynamic drag and the weight of the wing. 


Optimization techniques were applied in references 4 to 6 to minimize the mass 
of wing structures represented by finite elements. Reference 4 used the interior 
penalty function method of optimization with strength, frequency, and flutter as 
constraints. Reference 5 describes a computer program for the design of minimvim mass 
wing structures under flutter, strength, and minimum gage constraints. In refer- 
ence 6, a system of integrated computer programs is used for sizing structural member 
cross sections to satisfy strength and flutter design requirements for several vari- 
ants of an arrow wing supersonic cruise vehicle. The resulting structural weights 
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are used to determine the structural efficiency of the planform geometry, structural 
layout, type of construction, and type of material including composites. Several 
authors (refs. 7 and 8) have also considered the optimum design of structures with 
the material properties and loads treated as probabilistic quantities. However, 
probabilistic approaches to the design of airplane structures appear to be applied 
rarely. 

The present work presents a procedure for the automated optimum design of air- 
plane wing structures which limits the stresses developed due to landing and gust 
loads while considering static strength and flutter constraints. The dynamic loads 
resulting from the landing impact of airplanes contribute significantly to structural 
fatigue. In addition, the associated induced vibrations contribute to crew and pas- 
senger discomfort. Ground-induced dynamic loads and vibrations and the resulting 
problems will be magnified for supersonic cruise airplanes because of the increased 
structural flexibility of the slender body, thin-wing configuration and the higher 
take-off and landing speeds. Similarly, the dynamic stresses developed in an air- 
plane structure in flight due to atmospheric gusts are very important and often con- 
trol the design of wings of large airplanes. Since the contact conditions during 
landing ( such as sink speed, horizontal speed at the instant before contact, and 
rolling velocity) and the vertical velocity due to gusts are random in nature, exact 
values of the stresses induced in the airplane wing structure due to landing and 
gusts cannot be specified. Hence, the design optimization problem is solved twice, 
first by using a deterministic approach and then by using a probabilistic approach 
for the computation of the landing and gust stresses. The material properties are 
assumed to be deterministic in both the approaches. 

Two numerical examples are presented for demonstrating the feasibility of the 
method presented. First, the thickness and chord length of a uniform hollow double- 
wedge airfoil were determined with the help of graphical optimization techniques. In 
the second example, the thicknesses of cover skin and webs and the cross-sectional 
areas of flanges of a multiweb supersonic wing structure were found by using finite- 
element idealization and nonlinear programming techniques. The procedure described 
in this work is expected to be useful during preliminary design of an airplane. 


SYMBOLS 

a lift-curve slope 

c chord 

Cq reference chord 

D total drag 

friction drag 
Dp pressure drag 

f objective function 

F force transmitted by the landing gear 

Fg force between the tire and the ground 
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acceleration due to gravity 
jth constraint function 

amplitude of the frequency response function 
lift factor 

seitiispan of double-wedge airfoil (example 1) 

scale of turbulence 

lift distribution 

total number of constraints 

unsprung mass below shock strut 

mass distribution 

flight Mach number 

generalized mass in ith mode 

number of design variables 

number of segments of the airfoil (example 1) 
number of natural frequencies constrained 
atmospheric pressure 
generalized force in ith mode 

number of natural modes considered in the modal analysis 

standard deviation 

airfoil thickness 

atmospheric temperature 

vertical velocity of the gust 

flight velocity 

flutter speed 

forward velocity of flight at the time of gust encounter 
horizontal velocity at the instant before contact (during landing) 
sinking velocity of the airplane 
deflection of the unsprung mass 



w(x,y,t) deflection of the middle surface of the wing 

ith natural mode shape 

x,y coordinates in the middle plane of the wing 

Xj^ ith design variable 

x^/Y^ coordinates of the landing gear attachment point 

X design vector 

Yj^ ith random variable 

z vertical coordinate of wing 

tt/p constants in tire force~def lection relation 

variation of angle of attack of ith segment of the airfoil (example 1) 

a root angle of attack 

o 

Y ratio of specific heats 

6 wing tip deflection 

5^(t) ith generalized coordinate 

p density of air 

density of the material of the wing 
p^ solidity ratio 

a maximum gust stress 

g 

maximum landing stress 

a maximum steady- state stress 

s 

T time 

time at which the shock strut starts deflecting after initial contact of 
the tire with the ground 

duration of flight 

<])( t) random disturbance 

1 4>( t) Wagner function 

<^>{ w) power spectral density 

(|;(t) Kussner function; also random response of the system 


4 


<Hy 

1 z q 


function of random variables 


03^ ith natural frequency 

Superscripts: 

* value at the optimum point 

lower bound 
u upper bound 

Subscripts: 
b bending 

t torsional 

A bar over a symbol denotes mean value. A dot over a symbol denotes partial 
differential with respect to time. 


FORMULATION OF OPTIMIZATION PROBLEM AND SOLUTION PROCEDURE 

Loads Considered 

The following load conditions are considered in formulating the problem: 

1. The condition that exists during a longitudinal pull-up maneuver is consid- 
ered the static load condition. The root angle of attack and the deformed shape of 
the wing, which provide the specified gross lift, are determined. Then the stress 
distribution and, hence, the maximum principal stress at the wing root corresponding 
to the static displacement state are then computed. 

2. The landing load condition consists of the dynamic loads developed in the 
flexible wing structure during landing. Representative values of vertical (sink) 
velocity, forward velocity, and landing gear characteristics of the airplane are 
assumed. The interaction between the landing gear and the deformation of the air- 
plane structure is considered in finding the landing load history after initial 
touchdown. From this computed load history, the stress variation and, hence, the 
maximxim principal stress at the wing root during landing are determined. Since land- 
ing conditions such as sink velocity, forward velocity, and rolling speed are uncer- 
tain, they are treated as random variables in the probabilistic approach. 

3. The gust load condition consists of the loads developed in the flexible air- 
plane flying through a specified vertical gust. The vertical gust velocity distribu- 
tion is assumed to be cosine type with known maximum value. The time variation of 
the gust loads and, hence, the stress distribution in the aircraft wing is confuted 
by considering the gust magnitude to be constant in the spanwise direction. The 
random and continuous nature of the atmospheric turbulence is accounted for by assum- 
ing the gust velocity to be a stationary random process in the probabilistic 
approach. 
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Problem Formulation 


A general design optimization problem involves determining the design variable 
vector 


r 


= \ X* 


X 




( 1 ) 


which minimizes the objective function f(X) subject to the constraints 


g.(X) < 0 
D 


( j — 1 / • • • / m) 


(2) 


In this work, the objective function is the structural mass of the wing. In the 
first example, the total energy required to drive the airfoil through a specified 
flight condition is an alternative objective function. The following constraints are 
considered in the formulation of the problem. 

Natural frequency constraints .- The first p natural frequencies of the struc- 
ture 0 )^ are to lie between specified lower and upper bounds. Thus, 


g. (X) = (to5/w. ) - 1 < 0 


= (w./w’f) - 1 < 0 


1 1 ^ 


(i — 1, 2, .../ p) 


(i““ 1, 2, ..•/ p) 


J 


(3) 


Flutter constraint .- The flutter speed, Vp at a specified altitude is to be 
greater than a prescribed value Vp. Thus, 


®2p+1«> " ‘ " 


(4) 


Deflection constraint .-’ The deflection of the wing tip 6 under a specified 


u 


steady-state design condition is to be less than the maximum permissible value 6 . 
Thus, 


^2Ert.2«> - ‘ “ 


(5) 


Angle of attack constraint .- The root angle of attack of the wing a under a 

given steady-state design condition is to be less than a specified maximum value 

a . Thus , 
o 


^2p+3 



< 0 


( 6 ) 


6 


Steady-state stress constraint **- The maximum principal stress induced at the 

root of the wing a during a prescribed steady-state design condition is to be less 

than the maximum permissible value o • Thus, 

s 


^2p+4 



< 0 


(7) 


Dynamic stress constraints The maximum principal stress induced at the root of 
the wing during landing is to be less than the specified maximum value a ^ • 

Thus, ^ ^ 


^2p+5<^> ■ (V°t) - ' ‘ 

The maximum principal stress at the root of the wing due to specified gust conditions 
0 is to be less than the permissible value a • Thus, 

g g 

92p+6<*’ * ^ ° 

Equations (7) to (9) assume that the maximiim principal stress theory of failure 
is used instead of the more commonly used Von Mises theory. Further, the maximvim 
stress may be induced at a point other than the root of the wing. Hence, for a more 
accurate formulation of the stress constraints, the maximum stress obtained by con- 
sidering all the points of the structure is to be used in place of the maximum prin- 
cipal stress induced at the root of the wing in equations (7) to (9). 

Side constraints .- The design variables are to take values between speci- 

fied lower and upper limits. Thus, 

^2p+6+i^^^ = - 1 < 0 (i = 1, 2, n) (10) 


^2p+n+6+i^^^ = - 1 < 0 (i = 1, 2 n) 


Solution Procedure 

The nonlinear optimization problem stated in equations (1) and (2) is solved by 
the following techniques. 

1 . For the first example, which involves only two design variables (n = 2) , a 

graphical optimization procedure is used (ref. 9). 

2. For the second example, which involves more than two design variables, the 

Zoutendijk's method of feasible directions, coded in a program called 
CONMIN (ref. 10), is used for finding the optimtim solution. 
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ANALYSIS 


As stated earlier, the optimum design problem is solved by using a deterministic 
as well as a probabilistic approach for the coirputation of the stresses developed in 
the wing structure due to landing and gust loads. The uncertain parameters influenc- 
ing the landing and gust loads are assumed to be random variables with known mean 
values and standard deviations. In the deterministic approach, the mean values plus 
three corresponding standard deviations are used as the representative values for 
these parameters. These representative values will not be exceeded 99.97 percent of 
the time if the uncertain parameters follow the normal distribution. 


In the probabilistic approach, the mean values and standard deviations of the 
response quantities (stresses) are computed by using the mean values and standard 
deviations of the uncertain parameters, and the landing and gust stress constraints 
become 


u 




3S 


- 1 = 


u 


- 1 < 0 


( 12 ) 


and 


a 


u 

a 


9 


a 

g 


+ 3S 

a 

2 ^ 

u 


1 < 0 


(13) 


where the bar over a symbol denotes the mean and S the standard deviation of the 
landing stress and gust stress 1^ and follow the normal 

distribution, the constraints of equations (12) and (13) correspond to a confidence 
level of 99.97 percent. 


A brief outline of the landing and gust stress analyses according to determinis- 
tic and probabilistic approaches is given below. 


Deterministic Approach 

Landing analysis .- The stresses developed in the wing during landing are com- 
puted by considering the interaction between the landing gear and the flexible air- 
£)lane structure. The landing gear is assumed to have nonlinear characteristics 
typical of conventional gears, namely, velocity squared damping, polytropic air- 
compression springing, and exponential tire force-deflection characteristics. The 
derivation of the coupled nonlinear differential equations of motion that arise in 
the landing analysis and their numerical solution follow closely the procedure out- 
lined in reference 11. 
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According to the procedure described in reference 11, the displacement of the 
middle surface of the structure is expressed in terms of its natural modes of vibra- 
tion as 


w(x,y,T) = ^ W^^^(x,y) ^.(t) (14) 

i=0 ^ 

where W^^^(x,y) is the ith natural mode shape (i = 0 indicates the rigid body 
translation mode with W^^^(x,y) = 1) and is the ith generalized coordinate. 

Equation (14) shows that the first r natural modes are used in approximating the 
displacement of the airplane wing structure. If denotes the transverse 

displacement of the landing gear attachment point then from equation (14), 

one can write 


w^(t) = w(x^,y^,T) 


r 

Z 

i=0 


W 


(i) 


?.(T) 


(15) 


where denotes the value of W ^ (x,y) at equations of motion of 

the structure with the landing gear arrangement shown xn figure 1 can be derived 
in terms of' generalized coordinates as follows (ref. 11): 


M. 

X 




+ M. 0). (t) = Q. 

X X ^x X 


(x— 0, 1, 2, ..., r) 


(16) 


and 


- ff [l(x,y) - g m(x,y)] W^^Nx,y) dx dy 
Wing 
area 

(i = 0, 1, 2, ..., r) (17) 

where is the generalized mass of the structure in the ith mode (Mq is the 

total mass of the structure), (o. is the natural frequency of vibration in the ith 

mode (o) =0), is the generalized force in the ith mode, F is the force 

transmitted by the landing gear, L(x,y) is the lift distribution on the wing, 
m(x,y) is the mass distribution of the wing, and g is the acceleration due to 
gravi ty . 

Since the oleo-pneumatic shock struts are preloaded with air, they do not 

deflect until some time after the initial contact of the tire with the ground. 

d 

During the time 0 < t < t^, which can be called the first stage of impact, the shock 
strut is effectively rigid; thus the deflections of the unsprung mass and the landing 
gear attachment point can be taken to be same “ w^(i;)) for small values of 

0, the inclination of the shock strut. The motion of the complete system of landing 
gear and airplane is governed only by the force between the tire and the ground (Fg) . 

During the subsequent time t > which can be called the second stage of impact, 

the shock strut starts to deflect; the motions of the landing gear attachment point 
and the unsprung mass are different ('^j^(t) w^(t)). Hence, the equation of motion 

of the unsprung mass must be considered separately. Thus, the governing equations 
are different in the two stages. 
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For 0 < T < the governing equations are given by equation (16) with 


^ ^g ^u^U ** ^U^ 


( 18 ) 


and Fg can be expressed by an exponential variation of tire force with deflection: 

F = a(w (19) 

g u 


where a and p are constants and 

For X > T^/ both and w^j^(t) are treated as unknowns in equation (16). 

Since F is a nonlinear function of w^, which can be expressed in terms of 
using equation (14), equation (16) represents a set of r nonlinear differential 
equations. 


In this work, the equations of motion are solved by using only the first flex- 
ural mode in addition to the rigid body mode (i.e., r = 1). The initial conditions 
for the first phase of impact, 0 < t < x^, are assumed to be 


Cq(t=0) = 0 


( 20 ) 


where is the vertical (sink) velocity of the airplane at the time of initial 

contact with the ground. The total lift of the wing (b^otal^ assumed to be 


total 


= L(x,y) dx dy = W 

Wing 




( 21 ) 


area 


is the lift 


where W^otal half the gross weight of the airplane and 
factor. The initial conditions for the second phase of impact, t > T-, are the 
terminal conditions of the first phase. A step-by-step numerical solution procedure 
(ref. 11) is used to find the values of the generalized coordinates 
( t) at T = 0 , At, 2 At, . . • , where At is a small time increment. 


Once the values of known, the transverse displacement of 

the middle surface of the wing structure can be obtained as 


w(x,y,T) = ^0^^^ ^ ^ 


( 1 ) 


(x,y) l^(x) 


(t = 0, At, 2 At, 


( 22 ) 


With the known deflection w(x,y,T), the stresses induced in the wing structure at 
various instants of time can be found from the basic relations of structural mechan- 
ics. The maximum principal stress developed near the root of the wing (c!^) , during 
the total time interval considered for the analysis, can then be found without much 
difficulty. 
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Gust analysis To find the stresses developed in the wing due to a gust, the 
equations of motion of the aircraft are formulated with the vertical motion and wing 
bending modes as generalized coordinates. The gust is considered discrete and uni- 
form along the span of the wing. 

The transverse displacement of the middle surface of the wing can again be 
expressed by equation (14), and the equations of motion in terms of the generalized 
coordinates by equation (16). The generalized force in the ith mode is 

given by 


Q. = ff F(x,y,T) W (x,y) dx dy (23) 

Wing 
area 

where F(x,y,T) denotes the applied force distribution on the wing. By neglecting 
the effect of chordwise bending of the wing and using a strip type of analysis, equa- 
tion (23) can be expressed as (refs. 12 and 13) 


= / F(y,r) W (y) dy (24) 

Wing 
span 

with 


ape V 

F(y,T) - " ^ / w(y,x)[l - <))(x-x')] <ix' 

ape V ^ 

+ — I ^Kx-x') dx' 

^ Jo 


(25) 


where the time x is zero at the start of the gust penetration, a is the lift- 
curve slope, p is the density of air, c^ is the reference chord, Vg is the 
forward velocity of flight, u is the vertical velocity of the gust, [1 - <J)(x)] is 
the Wagner function, and 4 j(x) is the Kussner function. The lift-curve slope a 
may be chosen to correct approximately for aspect ratio and compressibility effects. 

With equations (23) to (25), equation (16) for i = 0 and 1 can be written, 
after some mathematical manipulation, as a set of two simultaneous integral- 
differential equations (similar to eqs. (16) and (19) of ref. 12). These equations 
of motion are solved for the generalized coordinates ?-|(x) at the dis- 

crete time stations x = 0, Ax, 2 Ax, ..., by using the numerical procedure outlined 
in reference 12 and by assuming that the airplane is in level flight before encoun- 
tering the gust. The variation of the vertical gust velocity with time is assiamed to 
be a cosine with known amplitude. Once the generalized coordinates 
5^(x) are known, the transverse displacement w(x,y,x), the stresses xnduced, and 
hence Cg, the maximum gust stress induced at the root of the wing, can be computed 
from the standard relations of structural mechanics. 
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Probabilistic Approach 


Landing analysis *- Results of measurements of contact (touchdown) conditions of 
several transport airplanes made during routine daylight operations in clear air are 
given in reference 14. It was found that the contact conditions such as sink speed/ 
horizontal speed at the instant before contact/ and rolling velocity/ although 
treated as constants in the deterministic analysis, are really random in nature. 
Since a function of several random variables is also random/ the maximum landing 
stress should be treated as a random variable. Let (\) denote a general function of 
the independent random variables (y-| / Y 2 f •••/ whose mean values 

y_) standard deviations , .../ \ are known; then the 


1 


approximate values of the mean and standard deviation of (\> can be determined as 
(ref. 15) 




(26) 


and 


S « 

4- 


ii 


Yi,y2..- 


j 


1/2 


(27) 


In this work, the sink speed (V^) and the horizontal speed at the instant before 
contact (V.) are taken as normally distributed and independent random parameters, so 
that ^ 




(28) 


and 




5V 






V 


9V, 


V /V^ ^ 
v' V 


V, 


1/2 


(29) 


The partial derivatives of the maximum landing stress appearing in equation (29) 
are evaluated by using a finite-difference scheme: 


SV 


a„(v + Av ,vj 

1'' V V 1'' 




V /Vn 
v' a 


AV 


(30) 


V 


and 


5o„ 


ev , 








AV, 


(31) 


12 


where AV^ and AV^ are small increments in the values of and and 

cTj^(v^ + AV^,Vj^), cTj^(v^/V^ + AVj^), and are computed according to the 

deterministic landing analysis procedure outlined earlier. 

Gust analysis *- The main drawback of the discrete gust approach is that it 
depends on a standard gust profile shape which cannot be derived on a logical basis 
by either theory or experiment and which cannot take into account the random and 
continuous nature of atmospheric turbulence. Hence, a probabilistic approach which 
permits description of the random atmospheric turbulence and the associated airplane 
response in analytic form, should be used. If the atmospheric turbulence is consid- 
ered a stationary random process, power spectral methods can be used for finding the 
root-mean- squared values of the stresses induced in the airplane wing structure. 

It is well-known (ref. 16) that if (})(t) represents a stationary random distur- 
bance such as the atmospheric vertical velocity and (^(t) denotes a system response, 
such as the induced stress, then the power spectrum of the response is given by 

3>^(0)) = |h(w) (32) 


spectrum of the disturbance, or input 
spectrum of the response, or output 

amplitude of the frequency response function, which is defined as the 
system response to a sinusoidal disturbance of frequency cd 

In equation (32) , the system is assumed to be linear and the turbulence one- 
dimensional (i.e., at any instant of time, the vertical gust velocity is constant 
along the wing span, or the scale of turbulence or the area under the correlation 
curve (ref. 17) is very large compared with the span of the wing). The spectrum of 
the response can be used to determine the root-mean-squared value of the response 
S , as 


where 


^. ( (o) 

^,(0)) 

4> 

|h(u) 


= r ^>,(w) dw (33) 

^ J-a> 'I' 


In the present work, the spectrum of the vertical gust velocity is taken as 
(Dryden’s model (ref. 18)) 


$.( 0 )) 


2L 

V 

g 


1 + (SL^co^/V^)! 

[l + (co^lVv^]J 


(34) 
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where S ^ is the root-mean-squared value of the gust velocity given by 


s? = I do) (35) 

and L is the scale of turbulence and o) is frequency. The integral of equa- 
tion (33) is evaluated by a nxamerical procedure involving the evaluation of the spec- 
tr\jm of the response for several values of the frequency. For each frequency o)/ the 
frequency response function H(a))/ and hence the spectrum of the response, is com- 
puted by finding the stress response at a point in the airplane wing due to a unit 
sinusoidal gust. Since this sinusoidal gust can be treated as discrete, the deter- 
ministic gust analysis procedure outlined earlier is used to find the transfer func- 
tion H(o))* Thus it becomes necessary to use the deterministic gust analysis a num- 
ber of times in one probabilistic gust analysis. 


RESULTS AND DISCUSSION 

The design of two example wings is considered for illustrating the optimization 
procedure. The details of analysis and the niimerical results are given below. 


Double-Wedge Airfoil (Example 1) 

The first example is the design of a hollow symmetric double-wedge airfoil shown 
in figure 2. The airfoil thickness (t) and chord length (c) are the design variables 
and the structural mass is the objective function. This design problem is solved by 
minimizing the energy required to drive the airfoil through a specified flight condi- 
tion as an alternative objective. The two objective functions are given by 




(36) 


and 


f^(X) = DUt^ 


(37) 


where is the density of the structural material, Ji is the semispan, p is 

the solidity ratio of the wing, D is the total drag, U is the flight velocity, 
and is the duration of the flight. The total drag can be expressed as (refs. 1 

and 3) 


D = Dp + Df 


(38) 


where 


and D£ represent the pressure and the friction drag with 


N 


= E 2yP M 


i=1 


a + a. + 
o 1 




N 


(39) 
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(40) 


D 


f 



cos ( a + a. ) 

O 1 



and 


T1(X) = 


0.259p M 

CO 

(1 + 0.13M'^) 


\log 


(1.51 X 10^)pjyi[T^(1 + 0.13M^) + 110.39](x + 


-2.584 
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T (1 + 0.13M ) 




In equations (39) to (41), p denotes the atmospheric pressure, y the ratio of 
specific heats, M the flight Mach nximber, N the number of segments of the air- 
foil, and the shear stress- 

The design data are given in table I. For the purpose of structural analysis, 
the semispan of the airfoil is divided into five equal segments, and the bending and 
torsional flexibility influence coefficients are used to find the steady-state 
stresses and the natural frequencies of the airfoil by using standard methods of 
structural analysis (ref. 1) . The first bending and the first torsional modes of the 
airfoil are used as generalized coordinates for the bending-torsion flutter analysis 
of the structure. The airfoil is assumed to fly at a Mach number of 2.0 during 
steady-state flight. Piston theory is used for determining the root angle of attack, 
the static loads producing tip deflection and steady-state stress, the flutter veloc- 
ity, and the aerodynamic drag of the airfoil. 

The optimum solutions are found graphically by plotting the contours of the 
objective functions (eqs. (36) and (37)) and the constraint boundaries. The contours 
of the objective functions in the design space are shown in figure 3. For the deter- 
ministic approach, the feasible design space and the optimum points are indicated in 
figure 4. If the structural mass of the airfoil is minimized, the optimum solution 
corresponds to point D.| in figure 4 with 



and f* = 1077.6 kg. This solution corresponds to an energy of 162.54 GJ. It can be 
seen that the upper bounds on the root angle of attack and the landing stress are 
active at the optimiim point. On the other hand, if the energy is minimized, the 
optimum solution is given by point D 2 with x = [0.246 m, 3.81 m] and f* = 80 GJ. 
This solution corresponds to a structural mass of 1422.0 kg. The lower bound on the 
torsional frequency and the upper bound on the landing stress are active at this 
optimum point. These results can be used to determine the mass penalty to be paid 
when minimizing energy or the energy penalty when minimizing structural mass. The 
point corresponds to t/c = 0.144, while the point D 2 corresponds to 
t/c = 0.065. A smaller value of t/c is required to reduce the energy, while a 
larger value is required to reduce the structural mass of the wing, as would be 
expected. 

For the probabilistic approach the contours of the objective functions, con- 
straint boundaries, and the feasible design space are shown in figure 5. The optimum 
solutions corresponding to minimum mass and minim\im energy are given by points and 


15 



p 2 f respectively. The point P-| corresponds to ^ = [0.360 m, 2.2 20 with a mass 

(minimum) of 1217.7 kg and an energy of 196.38 GJ. The point P^ represents the 
design vector X = [0.276 m, 4.305 m] which corresponds to a mass of 1810.15 kg 
and an energy (minimum) of 86.2 GJ. It can be seen that the points P-^ and P 2 
correspond to t/c equal to 0-162 and 0.064, respectively. As in the deterministic 
case, a smaller value of t/c is necessary to reduce the energy, while a larger 
value is necessary to reduce the structural mass of the wing. 

For the same value of the permissible landing and/or gust stress, the minimum 
value of the objective function is higher if probabilistic analysis is used instead 
of deterministic analysis- For example 1, the optimum points corresponding to mini- 
mum mass are (with f* = 1077.6 kg) for the deterministic analysis and P-j (with 
f* = 1217.7 kg) for the probabilistic analysis. Thus the minimum mass from the prob- 
abilistic analysis is 140.1 kg higher. 


Supersonic Transport Wing (Example 2) 

The design of the supersonic wing structure shown in figure 6 is the second 
example. The overall dimensions of this wing correspond approximately to those of 
the Concorde wing. The design data are shown in tables II and III. The wing is 
modeled with three types of finite elements. Constant strain triangular membrane 
elements represent the wing cover panel, quadrilateral shear web elements represent 
the ribs and spars, and pin- jointed bar elements represent the axial-load-carrying 
capacity of ribs and spars. The structural mass is the objective function. Six 
design variables are used in the optimization, the first four corresponding to the 
thickness of the skin in the four regions indicated in figure 6(a), the fifth one 
denoting the thickness of the rib and spar webs, and the sixth one representing the 
cross-sectional areas of the pin- jointed bars. 

A clamped boundary condition is specified along the root of the wing. Only the 
top half of the wing is idealized for the finite-element analysis by using equivalent 
plate elements. Since the load vector in the steady-state design condition depends 
on root angle of attack as well as on displacement state (i.e., nodal displacements), 
an iterative process was used to determine the correct values of the root angle of 
attack and the nodal displacements, which provide the specified gross lift (ref. 19). 
The stress state induced in the finite elements is determined from the known nodal 
displacements according to the stress-strain and the strain-displacement relations of 
linear elasticity. The maximum principal stress at the root of the wing a is then 
computed. 

The natural frequencies and mode shapes of the wing are found by setting up the 
eigenvalue problem and by reducing its order to one-third of that of the correspond- 
ing static problem. With Guyan's reduction technique (ref. 20), only the transverse 
displacement degrees of freedom are retained. 

The flutter problem is formulated by using the first four natural vibration 
modes as generalized coordinates and is solved by using a double iterative scheme 
(ref. 4). Piston theory (ref. 21) is used for computation of the root angle of 
attack, the aerodynamic loads producing tip deflection and steady-state stress, and 
the flutter velocity of the wing structure. The rigid body mode and first bending 
mode are used as generalized coordinates in the landing and gust stress analysis. 
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The number of degrees of freedom considered in the finite- element analysis and the 
number of modes considered in the dynamic analysis are kept low for computational 
convenience- 

A computer program was developed for implementing the present analysis (includ- 
ing the static, eigenvalue, and flutter analysis, and the deterministic or probabi- 
listic landing and gust stress analysis) in conjunction with the optimization pro- 
gram. This program was used for solving the two numerical examples. For example 2, 
the results obtained using the deterministic approach are shown in tables IV and V. 

It can be seen that the design variables X 2 r x^, and X 0 attained their lower 
bounds and x^ nearly reached its lower bound at the optimum point. Among the 
behavior constraints, the lower bound constraint on the second natural frequency was 
active and the flutter speed was near its lower bound. The reduction in structural 
mass obtained at various stages (iterations) of optimization is shown in figure 7. 

The reduction in the objective function is negligible after 12 iterations. 

In practice, the optimum values of design variables must be rounded off to the 
nearest practical values. In some cases the designer may be interested in altering 
the optimum design variables to satisfy certain new design requirements. When the 
various design variables are changed, the values of the objective function and the 
response quantities also change. Hence, a sensitivity analysis of the objective 
function and the behavior constraints with respect to changes in various design vari- 
ables was conducted at the optimiom point. The optimiam design was taken as the refer- 
ence design and the value of each of the design variables was changed by ±30 percent 
in steps of 10 percent. 

Figure 8 shows the results of this sensitivity analysis. The structural mass 

(fig. 8 (a)) is quite sensitive to variation of the design variable x^ . Figures 8 (d) 

and 8 (e) show that and a are also most sensitive to variations of x- - How- 

s g 1 

ever, (fig* 8 (b)) and Vp (fig. 8 (c)) are most sensitive to x^. For example, 

changing any one of the design variables x-^ , X 2 , X 3 , x^, or Xg within 

±30 percent of their optimum values does not violate the constraint on Vp,* but 
reducing x^ by more than 15 percent does violate this constraint. Figure 8 can 
also be used to identify the less sensitive design variables, so that the designer 
need not consider them while designing a similar system. 

The initial design vector used in the deterministic approach violated the proba- 
bilistic landing and gust stress constraints. Hence a different vector is used as 
the initial (feasible) design vector for the probabilistic optimization approach. 

All the lower and upper bounds on the response quantities and the design variables 
were the same. The results of optimization are shown in tables VI and VII. The 
minimum mass from the probabilistic analysis is 3645.4 kg, while the corresponding 
value from the deterministic analysis is 2858.4 kg. Thus the optimum value of the 
objective function is higher in the probabilistic approach. 

The design variables X 2 , X 3 , X 5 , and Xg assumed their lower bounds at the 

optimum point. Among the behavior constraints, the lower bound on the second natural 
frequency and the upper bound on the gust stress are active at the optimum point. 

Comparing the results obtained using the deterministic and probabilistic 
approaches leads to the following conclusions: 

1. The optimum value of the mass is higher if a probabilistic approach is used 
for the computation of the landing and gust stresses. 
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2. The set of active behavior constraints at the optimum point is different for 

the deterministic and probabilistic approaches- 

3. The design variables X 2 , / x^ , and x^ assxime their respective lower 

bounds for both approaches- 

These design variables may not reach their lower bounds if the maximum principal 
stresses used in the constraint equations (7) to (9) are determined by considering 
all points of the wing structure instead of only those at the root of the wing. 


CONCLUDING REMARKS 

A procedure is described for the automated optimum design of airplane wing 
structures with consideration of the stresses developed due to landing and gust loads 
along with the static strength and flutter constraints. Both deterministic and prob- 
abilistic approaches are used for finding the landing and gust stresses. The proce- 
dure is demonstrated by considering the design of two example wings: one based on 

simple beam-type analysis and the other based on finite-element analysis. 

The optimization results indicate that for the same value of the permissible 
landing or gust stress, the minimxim value of the objective function is higher if 
probabilistic analysis is used instead of deterministic analysis. Since the actual 
values of the stresses induced due to landing and gust loads cannot be predicted 
precisely, the probabilistic approach is expected to be more accurate and, hence, 
should be used in the optimum design of airplane wing structures. The procedure 
outlined is expected to be useful during preliminary design of an airplane structure. 

A sensitivity analysis was performed on the optimization results for the finite- 
element model. The results show that the structural mass, the maximum steady-state 
stress, and the maximum gust stress were most sensitive to variation of skin thick- 
ness near the wing root, while the second natural frequency and flutter speed were 
most sensitive to variation of skin thickness near the wing tip. These results are 
expected to be useful when the designer is interested in altering the optimum design 
variables to satisfy certain new design requirements. 


Langley Research Center 

National Aeronautics and Space Administration 
Hampton, VA 23665 
May 19, 1982 
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TABLE I.- DESIGN DATA FOR EXAMPLE 1 


Properties of wing: 

Semispan, • • • 9 

Solidity ratio of cross section 0.075 

Material (titanium) properties: 

Young' s modulus, GPa 110 

Shear modulus, GPa 44 

Density, kg/m 4500 

Additional masses: 

Mass of fuel (uniformly distributed in wing), kg 3000 

Mass of fuselage, kg 4000 

Mass of engine (concentrated 5.4 m from wing root), kg 2000 


Static load condition (pull-up maneuver): 

Altitude, m 10 668 

Air pressure. Pa 23 844.4 

Speed of sound, m/sec 296.7 

Air density, kg/m^ 0.3799351 

Flight velocity, m/sec 593.4 

Duration of flight, sec 5400 

Pull-up acceleration 2g 


Gust load condition: 

Discrete gust (cosine shape); 

Maximum vertical velocity, m/sec 3 

Length, chords 20 

Random gust: 

Standard deviation of vertical velocity, m/sec 1 

Scale of turbulence, m 762 

Forward velocity of flight, m/sec 154.8 

Altitude, m 7620 

Air density, kg/m^ 0.5494 


Landing load condition (fig. 1): 

Unsprung mass below shock strut, kg 

Landing gear distance from root of wing, m . 

Lift factor 

Duration of transient analysis, sec 

Tire force-deflection relation: 

a 

p 

Inclination of shock strut, rad 

For deterministic analysis: 

Sink velocity, V^, m/sec 

Horizontal velocity at contact, V^, m/sec 
For probabilistic analysis: 

V^, m/sec 



320 

1.8 

1.01 

0.3 

1.25 X 10^ 

1.22 

0 

3 

45 

.... 2.3 

.... 0.1V 

V 



m/sec 


34.6 

o.i5j 
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TABLE I • - Concluded 

Bounds on design variables and other constraints: 

Lower bound Upper bound 


=t, m 0*18 0 .45 

= c, m 1.5 4.5 

0)^ = co^, Hz 2.0 4.0 

(i) = / Hz 10*0 30*0 

6/ iti 1.2 

a / rad 0.08 

Vp, m/sec 890.0 

a , MPa 100 

Cw MPa 200 

MPa 800 

s 
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TABLE II.- DESIGN DATA FOR EXAMPLE 2 


Material (aluminum alloy) properties of wing: 


Young's modulus, GPa 68.948 

Poisson ' s ratio 0.3 

Density, kg/m^ 2768 

Additional masses: 

Maximum take-off mass, kg 174 632.9 

Mass of all engines, kg 5669.9 

Fuselage and payload, kg 32 658.6 

Fuel, kg 41 957.99 

Static load condition (pull-up maneuver) : 

Altitude, m 7620 

Air pressure. Pa 37 649.7 

Speed of sound, m/sec 309.7 

Air density, kg/m^ 0.5498382 

Flight velocity, m/sec 526.49 

Pull-up acceleration 3.75g 


Gust load condition: 

Discrete gust (cosine shape) : 


Maximum vertical velocity, m/sec 1.5240 

Length, chords 10 

Random gust; 

Standard deviation of vertical velocity, m/sec 0.5080 

Scale of turbulence, m 762 

Forward velocity of flight, m/sec 154.8 

Altitude, m 7620 

Air density, kg/m^ 0.5498382 


Landing load condition; 

Unsprung mass below shock strut, kg 
Landing gear location: 



™ 

Lift factor 

Duration of transient analysis, sec 
Tire force-deflection relation: 
a 


p 

Inclination of strut, rad 

For deterministic analysis: 

Sink velocity, V^, m/sec .... 
Horizontal velocity at contact 
For probabilistic analysis: 

V , m/sec 

"v 

_ V 

V^, m/sec 


m/sec 


Number of natural frequencies constrained 


2264.4 

19.3548 

2.7737 

1.01 

0.6 

16.8998 

1.21 

0 

1.1430 

45 

0.8792 

0.1V 

V 

34.6 

o.i5j 

2 
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TABLE III.- THICKNESS AND FUEL DISTRIBUTION 


FOR EXAMPLE 2 


z-coordinate 

Mass of 

of top node. 

fuel. 

Iran 

kg 

45.72 

1961.82 

35.56 

1519.56 

60.96 

1698.73 

96.01 

2694.38 

55.88 

1555.85 

83.36 

2036.66 

146.30 

3447.36 

111.25 

2599.13 

76.20 

1791.72 

121.92 

1258.74 

177.80 

1834.81 

205.74 

2122.85 

151.13 

1555.85 

96.52 

984.31 

167.64 

1950.48 

231.14 

2698.92 

299.01 

3946.32 

243.84 

2855.41 

180.34 

2088.83 

116.84 

1356.26 

223.52 

0 

297.18 

0 

373.38 

0 

411.48 

0 

320.04 

0 

228.60 

0 

137.16 

0 


Node numbers are defined in figure 6(b) 











TABLE IV.- DETERMINISTIC OPTIMIZATION RESULTS FOR EXAMPLE 2 



Initial 

design 

Bounds 

Optimiim 

design 

Lower 

Upper 

Design variables: 

, mm 

, mm 

x^, mm 

Xc / mm 

^ 2 

x^ , mm 

Structural mass of wing 

(objective function) / kg 

3.81 

3.81 

3.81 

3.81 

6.35 

161.3 

7323.3 

1.016 
1.016 
1.016 
1 .016 
1.016 
25.8 

12.7 

12.7 

12.7 

12.7 

12.7 

322.6 

2.7211 

^1.0160 

1.0244 

1.4917 

^1.0160 

^25.8 

2858.4 

Total number of function evaluations required 

116 

Computer time required,^ CPU sec 

3060 


^ower bound value. 

^Control Data CYBER 175 Computer. 
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TABLE V.- BEHAVIOR CONSTRAINTS FOR EXAMPLE 2 FROM DETERMINISTIC ANALYSIS 


Behavior constraint 

Initial 

Bounds 

Optimum 

design 

design 

Lower 

Upper 

Wing tip deflection, 

6/ m 

1.0158 


1.2700 

1.2147 

Root angle of attack, 

a , rad (deg) 

0.1539 


0.2088 

0.1601 

o 

First natural frequency, 

0)^ , Hz 

(8.8178) 

0.9615 

0.4775 

(11.9633) 

1.5915 

(9.1730) 

0.5777 

Second natural frequency, 

W2# Hz 

2.6666 

1.5915 

3.9789 

^1.5976 

Flutter speed, Vp, 

m/sec 

1944.7759 

526.5024 


601.7922 

Steady state stress, 

a , MPa 

s ' 

275.4886 


344.740 

®338.1 141 

Landing stress, cr., 

MPa i 

^413.1502 


413.688 

^409.1581 

Gust stress, c / MPa 

149.6999 


206.844 

^206.8440 

g 1 






^Active constraint* 
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TABLE VI.- PROBABILISTIC OPTIMIZATION RESULTS FOR EXAMPLE 2 



Initial 

design 

Bounds 

Optimum 

design 

Lower 

Upper 

Design variables: 

mm 

mm 

Xc # mm 

^ 2 

Xg, mm^ 

Structural mass of wing 

(objective function) / kg 

5.08 

3.81 

3.81 

3.81 

6.35 

161.3 

7994.5 

1.016 

1.016 

1.016 

1.016 

1.016 

25.8 

12.7 

12.7 

12.7 

12.7 

12.7 

322.6 

4.1783 

^1.016 

^1.016 

1.5519 

^1.016 

®25.8 

3645.4 

Total number of function evaluations required 

82 

Computer time required,^ CPU sec 

1648 


^ower bound value. 

^Control Data CYBER 175 Computer. 
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TABLE VII.- VALUES OF BEHAVIOR CONSTRAINTS FOR EXAMPLE 2 


FROM PROBABILISTIC ANALYSIS 


Behavior constraint 

Initial 

Bounds 

Optimum 

design 

Lower 

Upper 

design 

Wing tip deflection, 

6, m 

0.9806 


1.2700 

1.2236 

Root angle of attack, 

a , rad (deg) 

0.1537 


0.2088 

0.1593 

o 

(8,8063) 


(11.9633) 

(9.1272) 

First natural frequency, 

0)^ , Hz 

0.9888 

0.4775 

1.5915 

0.5917 

Second natural frequency, 

Hz 

2.7158 

1.5915 

3.9789 

^1.6325 

Flutter speed, Vp, 

m/sec 

1949.0350 

526.5024 


615.6526 

Steady- state stress, 

a / MPa 

s 

224.6464 


344.740 

238.5946 

Landing stress, o , 

MPa 

407.9998 

i 

413.688 

398.2367 

Gust stress, a , MPa 

159.8697 


206.844 

^206.4717 

g 






^Active constraint . 
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Figure 3.- Contours of objective functions/ structural mass, and energy 

for example 1 . 






Figure 5.- Graphical optimization of example 1 with probabilistic analysis 


used for a. 


and 














Figure 6.- Supersonic transport wing (example 2). 
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® Engine mass of 1417.5 kg at each of four locations 
Solid lines indicate spars and ribs 
Dotted lines indicate cover plate finite elements 

(b) Planform nodes and geometry. 

Figure 6.- Concluded. 
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0 2 4 6 8 10 12 14 16 

Iteration number 


Figure 7.- Convergence of objective function for example 2 with 

deterministic analysis. 




Change in x^, percent 


(a) Structural mass. 


Figure 8.- Sensitivity of various parameters to design variables for example 2 
with deterministic analysis used for and a * 



1.7027 



Change in percent 

(b) Second natural frequency • 


Figure 8.- Continued 



Change in x^, percent 


(c) Flutter velocity 


Figure 8.- Continued. 





a , MPa 
s 



Change in x^, percent 


(d) Maximum steady-state stress 
Figure 8.- Continued. 
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Change in x^, percent 


(e) Maximum gust stress 


a • 
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Figure 8.- Continued. 
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